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Abstract 

In view of the interpretation of data collected by large deep underground exper- 
iments in terms of primary Cosmic Ray physics, this work is focused on the study 



(N 
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OO \ of the production of of TeV muons in Extensive Air Showers. The review tries to 

point out those features of hadronic interactions that mostly affect the production 
of the high energy muons. A few different Monte Carlo codes are compared, with a 
\ particular attention to those based on the Regge-Gribov framework. The possibility 

of performing experimental tests of the proposed models is also briefly discussed. 
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1 Introduction 

Cosmic rays of the highest energies can be observed only through their interaction with the 
Earth's atmosphere: the reliability of our interpretation of the features of such secondary 
particles, and of their relation to the characteristics of the primary particle, is necessarily 
related to quality of our understanding of hadron-hadron, hadron-nucleus and nucleus- 
nucleus interactions. This aspect is particularly stimulating for high energy physicists, 
since there is not yet an exact way to calculate the properties of the bulk of hadronic 
interactions. Also, from the experimental point of view, the productions of secondary 
cosmic rays at very high energies occurs in kinematic regions, or energy ranges, that have 
not been explored in accelerator experiments, and that will hardly be accessed even at 
the hadron colliders of the next generation. Therefore, Extensive Air Showers of energy 
> 10 15 still represent an almost unique chance to test our theoretical achievements in very 
high energy nuclear physics, not to speak of the region above 10 17 eV, which is, at least 
now, completely out of reach of the present accelerator technology. 

In the last years, new experiments, located on the earth surface or underground, have 
collected valuable data concerning different Extensive Air Shower (EAS) components. 
The question of the interpretation of these data, the choice of the best possible simulation 
tool, and in general the evaluation of the systematics associated to such simulations, has 
become crucial. This remains a fundamental point also for the future activities at the 
extreme high energies 
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A decisive contribution to this discussion has come from the work carried on in 
the Karlsruhe group 0, concerning the maintenance and distribution of the CORSIKA 
code||, which now can be interfaced to different hadronic packages. A detailed com- 
parison of the results from these models has been presented in ||, and new progresses 
have been also presented at this conference||. There, the emphasis was given to the data 
collected by surface arrays, mostly the e.m. component and the low energy (Gev) muons. 

The purpose of this work is to present a further complement to the discussion started 
in ref.[[|], focusing the review on the Monte Carlo predictions which are specific for the 
production of very high energy muons (in the TeV range), as those collected by the large 
underground experiments. This is topic is of particular relevance in the indirect analysis 
of EAS, since TeV muons come from the decay of mesons produced in the very early stages 
of the cascade, and are more directly related to the properties of high energy interaction 
with respect to the low energy muon component. 

After the pioneristic work of Utah's experiment || and the results from NUSEX[[7] and 
FREJUS||, recent experimental papers on the subject have been published by experi- 
ments at Gran Sasso, like MACRO 0, [TO], and also by Soudan2|TT| and KGF|T2|. The 



MACRO data are important for the large collection area and the statistics, apart from 



the interest in the analysis method proposed in fllP] . The relevance of KGF results is 
mainly the high energy selection because of the large rock overburden. Both experiment 
based their analyses on Monte Carlo models based on parameterizations experimental 
results at colliders. A critical discussion of these models is now considered mandatory. 
Some of these experiments, like MACRO or LVD|]13]], take also data in coincidence [|T5| 



with the EAS-TOP surface array For this reason, the calculation of the TeV muon 
component has to be studied also in correlation with other EAS components. This will 
be the subject of a much more complete review [TO], of which the present work is just a 
partial anticipation. 

In the next section the features of hadronic interactions which are relevant for muon 
production are reviewed. Then, in Section 3, the Monte Carlo generators under study are 
briefly described. Some results about the simulation of primary interactions are given in 
Section 4, while the comparison of full shower calculations are presented in Section 5. A 
general discussion is attempted in the Conclusions. 



2 High energy muons as a probe of hadronic interac- 
tions 

The high energy muons in EAS come from the decay of pions and kaons produced in 
the early stages of shower development. The interest is therefore in the cross sections, 
multiplicities and their fluctuations, for meson production in all kinds of possible hadron 
and nucleus interactions in the relevant energy range. A key aspect is the competition 
of decay and interaction of these mesons along their path in a medium with varying 
density profile, such as the atmosphere. Therefore, while the inelastic cross section of 
primaries determines the height of first interaction, the inelastic cross sections of mesons 
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determine the interaction length. The decay length is a function of meson energy, whose 
distribution is governed by the Xp distribution as resulting from meson, nucleon or nucleus 
collisions with air nuclei. In particular, the high Xp values ("fragmentation region"), 
poorly studied in accelerator experiments, are the dominant ones. Diffractive and non 
diffractive interactions are both to be considered. 

The P± distribution affects the transverse structure of muons in the showers, since 
the kink in the meson decay can be considered negligible. At TeV energies, geomagnetic 
deflection is a minor correction. A realistic calculation must take into account the corre- 
lations existing between P±, energy and rapidity (or Xp). As already discussed in HT7f , 



while a simultaneous change in the interaction lengths of mesons and nucleons produces 
only second order effects in the eventual muon yield, variations in P± are directly reflected 
in the calculation of the lateral distribution. 

Of course, the lateral distribution is always the result of a convolution between longi- 
tudinal and transverse properties of the interaction. Calling r the displacement of a muon 
with respect to the shower axis, as measured at a slant distance H pro d from the relevant 
interaction, then: 

p 

r ~ — — H prod (1) 

In this simplified description the transverse momentum in the parent decay is neglected. 
The previous expression can be re-expressed in a more instructive way, considering that 
at high energy the longitudinal cm. variable Xp is approximately equal to the laboratory 
energy fraction X[ a b in the forward region, up to terms of the order of (m±/ Eq) 2 ) : 

oc 3^ (log a'^Mr + const.) (2) 

where a simple exponential atmosphere profile has been used in the last expression. 

Of course, in the case of underground experiments, muon propagation in the rock is 
also a fundamental aspect of calculation, but it is a factorizable ingredient, not to be 
discussed in this work. A discussion of this topic can be found in It can just be 

mentioned that, even for large depths, although the scattering in the rock is quantitatively 
important, the simulation results show that the amount of this deviations it is not yet 
enough to obscure the influence of P± in the lateral distribution of the surviving muon 
component. 

The energy filter provided by the rock, imposes a threshold on the muon energy, and 
therefore also on its parents' energy. This introduces some selection effects. Considering 
for example the case of Gran Sasso laboratory. There, the minimum rock overburden is 
about 3100 hg/cm 2 (in the direction of Campo Imperatore, where EAS-TOP is located). 
That depth corresponds to an energy threshold for muons of about E 1 ^ ~ 1.3 TeV. For 
an integral spectral index of primaries 7=1.7, the corresponding average energy of parent 
pions is 

nj _|_ 9 1 — r ( 7+1 ) 
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where r^=m 2 /m\. When considering the first primary interaction of a nucleon, at energy 
E , a threshold in muon energy is thus equivalent to a threshold in Xp ~ Xi a b = jjF, which 
for "low" Eq values (a few TeV) can be close to 1, i.e. in the very forward region. In that 
cases, in practice only muons from the first generation in the cascade are detected, and 
the primary interaction is really investigated. At fixed threshold, for increasing primary 
energy, the average X iab of parents at first decreases, and then exhibits a smooth rise, 
when more and more cascade generation can provide mesons above threshold. This is 
shown in Fig.|l| for parent pions from proton and Fe primaries, in the case of 3400 hg/cm 2 , 
at 30° of zenith angle. 

Since below 10 TeV, the composition of primaries contributing to TeV muons is 
strongly dominated by protons, this selection effect is being exploited in the coincidence 
experiment of EAS-TOP and MACRO to investigate the differences between interaction 



models in the fragmentation region] IS 



In this limited discussion the interesting question of prompt muons from the decay of 
heavy flavored mesons will be neglected. However, it is recognized that this production 
is indeed a very small fraction of TeV muons surviving underground (see for instance ref. 
26fl), and there is debate on the real possibility of a positive detection of this component. 



3 Monte Carlo generators 



After the first fundamental work of Elbert POU , an important paper concerning the simula- 
tion of muons detected deep underground is that of ref. plj , where parameterizations as a 



function of energy per nucleon, zenith angle and rock depth were presented, following the 

The calculations were performed using, as interaction 
algorithm" [p2f . This algorithm is indeed simple, perfectly 



method already pointed out inpO 
model,the "splitting Hillas' 
preserving Feynman scaling and does not provide any information about the transverse 
structure of the showers. More realistic predictions were obtained with the development 
of the HEMAS code[^3|], which was explicitly designed to provide a dedicated tool for the 
high energy muon component (above 500 GeV), in view of the experimental activity at 
Gran Sasso laboratory. The hadronic interaction model of HEMAS is based on the pa- 
rameterizations of collider data, mostly UA5[25|, following a picture in which clusters of 



particles were formed in the nucleon-nucleon collisions, eventually decaying into ordinary 
mesons. Effects due to the nuclear nature of the target, as measured in heavy ion physics, 
are used to correct multiplicity, transverse momentum, etc. 

The most noticeable differences coming from these collider parameterizations concern 
two aspects: 

1. the larger fluctuations of muon multiplicity, in connection to the negative binomial 
multiplicity distribution of produced mesons 

2. the transverse distribution of muons, following the power law distribution related 
to that of P± and its correlation with energy and multiplicity. 
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Figure 1: Monte Carlo calculation of Average Xi a b of parent pions of muons detected after 
3400 hg/cm 2 of standard rock, at 30° of zenith angle, for proton (black circles) and iron 
nuclei (open squares) primaries. 
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The first analyses by MACRO [0, 123] could establish the much better description of data 



provided by HEMAS with respect to the previous parameterizations of |^TJ. It must be 



quoted that also the recent analysis of KGF made use of similar collider parameterizations. 

In any case, this type of codes have still a weak point. The application to cosmic 
rays necessarily requires some extrapolation, in particular towards the high Xp values, 
and this procedure has always some degree of arbitrariness. Also, in the construction of 
parameterizations, important correlations among variables might be lost. These are some 
of the reason that have driven people to look for "theoretically inspired" models, with 
a limited number of parameters and a phenomenological framework as a guidance. A 
rather remarkable success in hadronic physics seems to be achieved by the Dual Parton 
Model [^| (DPM). This approach makes use of concepts derived from the mathematical 
requirements of scattering theory, as unitarity and analyticity together with color flow and 
parton idea and topological expansion of QCD. Such concepts can be formulated in the 



framework of relativist quantum field theory, according to the Regge-Gribov approach [^5 
The extension from nucleon-nucleon to nucleon-nucleus and nucleus-nucleus is achieved 
according to the Glauber multiple scattering formalism ||29||. A last fundamental ingredient 



of all these models is the "hadronization" algorithm, where parton strings produced in 
the scattering process are converted into ordinary particles. Some authors make use of 



the hadronization Lund models 30 



In practice, many different variations exist of these basic ideas, and also many different 
practical implementation have been developed in different numerical codes. Some of these 
have been successfully interfaced to the general purpose shower code CORSIKA, as already 
mentioned in the introduction. A detailed description of these models can be found only 
in the original literature, however the most important features are the followingr]. 

1. HDPM||. This is the original interaction package contained in CORSIKA, where 
parameterizations of the main results of DPM are contained. It will not be consid- 
ered in the following. 



SIBYLL||32||. It was explicitly developed for cosmic ray application. It contains DPM 
ideas, but the present version is deeply based on the minijet production as calculable 
from QCD, providing all the rising part as a function of energy of the nucleon- 
nucleon cross section. The soft part is constant. A Lund code for hadronization is 
used. 



3. DPMJET|32[]. It is based on the "Two component" DPM, i.e. both soft processes 
and hard QCD scattering are unitarized together. The soft part assumes the "super- 
critical pomeron" so that the soft contribution of the cross section is not constant 
but rises as s 08 . Multi-Pomeron exchange is possible. A specific point of merit 
of this code is the inclusion of dedicated algorithms to treat nuclear effects. Also 
DPMJET makes use of a Lund hadronization code. It must be mentioned that 
this interaction model has been also chosen for a new version of HEMAS [|35|1 , called 



1 CORSIKA makes use of the GHEISHA|3l| parameterized model below 80 GeV incident hadron 
energy. This is irrelevant for TeV muons. 
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HEMAS-DPM, again specifically devoted to high energy muons. This is presently 
in use by the MACRO experiment, and preliminary results are presented at this 
conference [pi | . From the point of view of results, no significant differences with 



respect to those achieved inside CORSIKA are expected, although the comparison 
of results from the same interaction model from two different shower codes would 
be a technically interesting issue. This will be discussed in [|i~6f. 



4. QGSJET|36|]. It is based on the "Quark Gluon String" model[[T?|]. It is very sim- 
ilar to DPMJET as far the pomeron parameters are concerned. It also includes 
the contribution of perturbative QCD. At present there are important differences 
with respect to DPMJET in the way in which soft and hard scattering are merged 
together. 

5. VENUS [^]. Originally developed for accelerator physics. It contains DPM with the 
supercritical pomeron, and the addition of some specific features concerning color 
exchange processes. It introduces secondary interaction of strings before hadroniza- 
tion. The version at present in CORSIKA does not yet include the contribution of 
minijets. 

In the following results concerning high energy muons from some of these choices will 
be compared. For some specific issues the with the original HEMAS code will be also 
shown. 



4 Particle Production in p— Air single interaction 

The first approach in the comparison of models for shower calculation is the analysis of 
the features of primary interaction. In principle, one has to study the primary interaction 
in both nucleon-Air and nucleus-Air collisions. However, for simplicity, the nucleon- 
Air interaction can be considered as the basic representative ingredient of the models. 
The extension to nucleus-nucleus interaction according to the Glauber model could be, in 
principle, a common part to all codes. However, a check of how this is realized in practice is 
mandatory and will be addressed in future. For the comprehension of shower development, 
the study of meson-Air interaction would be also important. An exhaustive examination 
would require a ponderous dissertation, that goes beyond the present possibility. The 
discussion will be limited to the inelastic cross section, the X[ ab and the P± distributions. 

4.1 The inelastic cross section 

In Fig. [|, the p-Air inelastic cross section as a function of laboratory energy is shown for 
some of the considered models, in the relevant energy range accessed with muons detected 
deep underground. 

There are significative differences in the knee region (VENUS seems to depart more 
from the other ones), but also at yfs values where good p-p data are available. Here 
the SIBYLL cross section seems to exhibit the major deviation. Presumably, this can be 
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Figure 2: Inelastic p-Air cross section as a function of laboratory energy for some of the 
considered Monte Carlo models. 
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attributed to a non correct inclusion of diffractive events in the Glauber calculation [|3S 
and in the next version of SIBYLL this problem will be corrected. These differences in 
cross section will be reflected by different average interaction heights and some change in 
the muon production height. 

4.2 Inclusive distributions and the Spectrum— weighted moments 

As a significant example of possible difference between codes, Fig. § shows the Xi ab 
distribution of charged pions produced in p-Air collisions at 200 TeV in the laboratory 
frame, for 3 different codes: DPMJET (continuous line), QGSJET (dashed) and VENUS 
(dotted). The difference of DPMJET with respect the other codes at large X is due to a 
particular mechanism of Lund fragmentation (JETSET 7.3) of the fast valence "diquark"Q 
in the projectile, in the case of positive pions ("popcorn" effect ^0f). No effect of this kind 
is instead possible for kaon production, where strange quark are not presented in the 
valence quarks, or diquarks, of the incoming nucleon. This mechanism is not included in 
the other MC codes. This feature has indeed relevance in the yield of high energy muon. 
It is important to notice that, apart from the region near X iab =0 and X lab =l, these 
Regge-Gribov model exhibit a substantial Feynman scaling. Therefore these features are 
preserved in a wide range of energies. 

In the current literature, some emphasis has been given to the "spectrum-weighted" 
moments, defined, for instance in case of pion production, as 

1 A l\Tp+Air—+n+X 

(x lab y<- 1 — — dx lab (4) 

where 7 is spectral index of primary nucleons. Such a function give a measurement of 
the inclusive particle yield in cosmic ray showers after the integration over the energy 
spectrum [[IT] . In principle, the high energy muon yield is proportional to those values. 



However, the comparison of Z-moment from different codes, as that shown in Fig. [|, can- 
not give completely meaningful information when the inelastic cross sections are different. 
In fact, a larger Z factor (as in the case of SIBYLL) will not give a larger muon yield if 
the production height is lower, so that available length for decay is smaller. 

It is also interesting to inspect the behaviour of the energy fraction Kb carried away 
by the leading baryon in p-Air collisions as a function of primary energy. This is shown 
in Fig.[5|. It can be seen how the HEMAS code exhibits a different behaviour from other 
representative models of the QGS and DPM kind, where elasticity smoothly decreases 
with energy. In this respect, it is worthwhile to remark what stated in ref . [|42 1 about 



the cluster models. Such approaches have the general tendency to produce events which 
become increasingly elastic as energy increase?]. In the same reference it was shown how 
QGSJET and DPM model represent interactions on nuclear targets in a way which is 

2 i.e. the projectile after the stripping of one valence quark 

3 It has been suggested that a less unbiased way to look at elasticity would be to consider K B _g, 
where anti-baryons are taken with the minus sign. 
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Figure 3: Xi a b distribution for charged pions produced in p-Air collisions at 200 TeV 
in the laboratory frame for DPMJET (continuous line), QGSJET (dashed) and VENUS 
(dotted). 
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Figure 4: Comparison of Z moments for charged pion production as a function of energy, 
for different codes. 
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consistent with low energy data. It can be important to notice how a model like the 
quoted Hillas' algorithm, in its original form, would give a constant Kb=0.5. 

A few remarks are necessary as far as the P± distribution is concerned. The phe- 
nomenological framework of the Regge-Gribov theories, by its nature, provides only pre- 
dictions for the longitudinal properties of the interaction. The transverse structure leading 
to the specific P± distribution is not constrained by the theory, but for the higher P± phe- 
nomena, where perturbative QCD can be used (of small relevance in the primary energy 
region addressed in general by the underground experiment). Once again, the model 
builders have to be guided mostly by experimental data, introducing a-priori determined 
functional forms with their additional parameters. Again, as a typical example, the P± 
distribution of charged pions produced in p-Air collisions at 200 TeV laboratory energy, 
as obtained in different models is shown in Fig.|6|. 

It can be pointed out how all the considered models produce P± distributions which, 
following the experimental results at colliders, can be fitted with a power law spectrum. 
It can be also noticed how HEMAS and DPMJET give very close results, predicting larger 
tails than the other models at relatively high values of P±. 



5 Shower Calculations 

Test simulation runs of EAS development have been performed, choosing typical param- 
eters for the experimental situation at Gran Sasso laboratory. According to the logic 
followed so far, results are presented here only for the case of proton primaries. Runs 
have been performed Runs were performed at different fixed energies, for fixed angles 
(30° zenith and 190° azimuth), setting the observation level for EAS at 2000 m a.s.L, 
and with the geomagnetic field specific for the site. The chosen direction corresponds to 
the typical one from which the underground halls see the EAS-TOP array. Propagation 
of TeV muons in rock has been calculated by means of the PROPMU code|33| for the 
corresponding thickness of 3400 hg/cm 2 of standard rock. In Table [l] the average yield 
of muons as measured at the surface (E M >1 TeV) is reported, together with that of the 
survived muons underground. 

A comparison of muon yields in showers from protons and primary nuclei is shown, 
for example at 2000 TeV total energy, in Table [|. The shape of the muon multiplicity 
distributions from the different models are very similar, as, for instance, in the example 
of Fig. |7], from p-Air collision at 2000 TeV. 

As far as the underground muon yield is concerned, it can be seen how these predictions 
are in general close one to the other. The most relevant percentage differences, for a given 
threshold on muon energy, are noticeable near the corresponding threshold energy for 
primaries. In this context, it can be important to remind the considerations expressed 
in Sect. 2, about the importance of high Xp region in meson production, where the 
distribution from different models depart more one from another. 

A summary concerning the quantities affecting the transverse structure of the high 
energy muon component is reported in Table |3[ These are the average depth of the first 
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Figure 5: Comparison of the energy fraction carried away by leading nucleons in p-Air 
collisions as a function of laboratory energy, for different models. 



13 




Figure 6: Comparison of P± distributions of charged pions produced in p-Air collisions 
at 200 TeV laboratory energy, for different models. 



14 



Table 1: Comparison of . 


average TeV muon yield at 2000 


m a.s.l. and of the no. of 


muons surviving underground, after 3400 hg/cm 2 , produced in proton induced showers 


at different energies. The 


same muon transport code has been applied in all runs. The 


relative statistical error on the reported figures is of the order of 5 


p-Air, 3 TeV 


Code 


< N^(E > 1 TeV) > 


< > survived undergr. 


HEMAS 


0.0027 


0.00049 


CORSIKA/DPMJET 


0.0035 


0.00059 


CORSIKA/QGSJET 


0.0018 


0.00015 


CORSIKA/VENUS 


0.0028 


0.00031 


CORSIKA/SIBYLL 


0.0042 


0.00056 


p-Air, 10 TeV 


Code 


< N^(E > 1 TeV) > 


< > survived undergr. 


HEMAS 


0.064 


0.015 


CORSIKA/DPMJET 


0.055 


0.016 


CORSIKA/QGSJET 


0.056 


0.015 


CORSIKA/VENUS 


0.057 


0.015 


CORSIKA/SIBYLL 


0.064 


0.018 


p-Air, 20 TeV 


Code 


< N„(E > 1 TeV) > 


< Nfj, > survived undergr. 


HEMAS 


0.182 


0.048 


CORSIKA/DPMJET 


0.142 


0.049 


CORSIKA/QGSJET 


0.153 


0.052 


CORSIKA/VENUS 


0.150 


0.049 


CORSIKA/SIBYLL 


0.153 


0.053 


p-Air, 200 TeV 


Code 


< N^E > 1 TeV) > 


< > survived undergr. 


HEMAS 


1.29 


0.47 


CORSIKA/DPMJET 


1.17 


0.49 


CORSIKA/QGSJET 


1.17 


0.48 


CORSIKA/VENUS 


1.28 


0.52 


CORSIKA/SIBYLL 


1.11 


0.44 
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Table 2: Comparison of average TeV muon yield at 2000 m a.s.l. and of the no. of muons 
surviving underground, after 3400 hg/cm 2 , produced in showers induced by different pri- 
mary nuclei at 2000 TeV total energy. 



p-Air, 2000 TeV 


Code 


< N^E > 1 TeV) > 


< > survived undergr. 


HEMAS 


6.44 


2.55 


CORSIKA/DPMJET 


6.90 


2.84 


CORSIKA/QGSJET 


6.79 


2.80 


CORSIKA/VENUS 


7.43 


3.09 


CORSIKA/SIBYLL 


5.99 


2.48 


He-Air, 2000 TeV 


Code 


< N^(E > 1 TeV) > 


< > survived undergr. 


HEMAS 


9.51 


3.75 


CORSIKA/DPMJET 


9.61 


3.94 


CORSIKA/QGSJET 


8.89 


3.66 


CORSIKA/VENUS 


10.26 


4.23 


CORSIKA/SIBYLL 


8.67 


3.59 


Fe-Air, 2000 TeV 


Code 


< N^E > 1 TeV) > 


< > survived undergr. 


HEMAS 


16.42 


5.61 


CORSIKA/DPMJET 


15.72 


5.86 


CORSIKA/QGSJET 


15.87 


5.74 


CORSIKA/VENUS 


16.92 


6.17 


CORSIKA/SIBYLL 


15.15 


5.57 
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Figure 7: Comparison of normalized distributions of underground muon multiplicity 
from p-Air collisions at 2000 TeV laboratory energy, 30o zenith angle, and 3400 hg/cm 2 
depth, for different models. 
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interaction Xfi rs t, < P± > for pions coming from the first interaction, the average slant 
production height of muons survived underground (the decay height of their parent 
mesonsfj), the average distance of muons from shower axis underground (< R >) and the 
average muon pair separation for multiple muon events (< D >). 

In the energy range around hundreds of TeV, to which most of data taken at Gran 
Sasso laboratory belong, the resulting differences in the average muon separation do not 
exceed 20%. These discrepancies seem to reduce at higher energy, while they appear much 
larger a few tens of TeV. DPM JET is probably the only model predicting a higher average 
separation than HEMAS. A precise analysis of the reasons leading to those differences 
in the models is complicated, however, it is important to notice that, in fact, HEMAS 
gives in general higher values of average P± with respect to the other models. The only 
exception is indeed DPMJET, which, as mentioned before, has a particular attention to 
the reproduction of nuclear effects affecting the transverse momentum, as measured in 
heavy ion experiments. On the other hand, the effect on this larger P± on the lateral 
distribution of muons is moderated in HEMAS by a deeper penetration of showers, giving 
in general a somewhat smaller average height of meson production. Similar features in 
the comparison of models are also obtained for nuclear projectiles. The data analysed by 



MACRO [[24], pTOH favour the HEMAS simulation, however a comparison of the same data 



with other models is still missing. 



6 Conclusions 

It has been shown that, as far as TeV muons are concerned, differences in the predicted 
yield are not enormous, although they exist. So far, the published analyses from different 
experiments have been mostly based on a given specific model. In the light of the present 
results, it is instead recommended that each group should consider more than one model 
for the interpretation in terms of spectrum and composition of primaries. In particular, it 
would be highly desirable that the set of considered simulation codes for event generation 
would be the same for every collaboration. This of course requires a considerable effort to 
the experimental groups, however it appears necessary in view of a sensible comparison of 
the different results, which at present are sometimes still controversial, despite the high 
quality reached in the detector and analysis techniques. This is valid not only for the 
underground experiments, but also for EAS array. In this respect, more complete results 



on code comparison will be provided in a next work [16]. The results presented at this 
conference by the MACRO experiment p4| , are an important example in the in the case 
of DPMJET, which provides a larger muon yield underground, suggested direction: the 
analysis of primary spectrum, originally based on HEMAS, can converge to significantly 
different values in the case of DPMJET, which provides a larger muon yield underground. 

A fundamental question is if it is possible an experimental discrimination of hadronic 
interaction models on the basis of cosmic ray data themselves. As far as high energy muons 



4 CORSIKA does not allow to access directly to the production height of parent mesons, which would 
be more interesting for our purposes 
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Table 3: Comparison of a few relevant quantities concerning the lateral distribution of 
underground muons at the depth of 3400 hg/cm 2 , from proton primaries at 20, 200 and 
2000 TeV, 30° zenith angle. 



p-Air, 20 TeV 


Code 


< X f irst > 


< Pj_ > 7T ± 


<H IM > 
(km) 


< R > 


< D > 




(g/cm 2 ) 


(GeV/c) 


(m) 


(m) 


HEMAS 


51.4 


0.40 


24.1 


7.9 


12.7 


CORSIKA/DPMJET 


44.4 


0.42 


25.6 


10.1 


13.9 


CORSIKA/QGSJET 


45.7 


0.39 


24.3 


7.3 


10.0 


CORSIKA/VENUS 


48.3 


0.35 


24.5 


7.4 


8.3 


CORSIKA/SIBYLL 


50.9 


0.37 


23.5 


7.2 


11.5 


p-Air, 200 TeV 


Code 


< A fi rs t > 


< Pj_ > 7T ± 


<H IX > 
(km) 


< R > 


< D > 




(g/cm 2 ) 


(GeV/c) 


(m) 


(m) 


HEMAS 


56.1 


0.44 


20.6 


5.3 


8.0 


CORSIKA/DPMJET 


53.9 


0.43 


21.7 


6.2 


8.8 


CORSIKA/QGSJET 


52.8 


0.41 


21.4 


5.5 


7.8 


CORSIKA/VENUS 


60.2 


0.36 


20.9 


5.3 


7.5 


CORSIKA/SIBYLL 


55.2 


0.41 


20.2 


5.2 


7.3 


p-Air, 2000 TeV 


Code 


< A f^st > 


< P ± > 7T ± 


<H IX > 


< R > 


< D > 




(g/cm 2 ) 


(GeV/c) 


(km) 


(m) 


(m) 


HEMAS 


63.0 


0.50 


16.3 


4.1 


6.0 


CORSIKA/DPMJET 


60.0 


0.42 


18.5 


4.9 


6.4 


CORSIKA/QGSJET 


63.1 


0.44 


17.7 


4.2 


5.6 


CORSIKA/VENUS 


66.7 


0.36 


16.8 


4.1 


5.3 


CORSIKA/SIBYLL 


60.3 


0.44 


17.0 


4.4 


5.6 
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are concerned, the answer is not yet certain. However, there are promising interesting 
measurements which cover the range of primary energy below few tens of TeV. The 



coincidence experiment Ne-N M at Gran Sasso[19] has been already mentioned. The same 



collaborations are now preparing the analysis of coincidences between Cherenkov signal 



and TeV muons underground ||44j| . This would allow the study of TeV muon yield as a 
function of energy in the threshold region. It is important to remark that, although these 
check can be done in a limited and relatively low energy, they are in any case important 
to understand features, like the shape of Xi a t, distribution which, as discussed before, are 
essentially preserved in a wide energy range. 

Last but not least, the quality of the different models can be debated just on phe- 
nomenological bases. Despite the common basic language, the presently available "the- 
oretically inspired" models do still contain many free parameters and rather different 
choices for the algorithmic solutions. Hopefully, if future data on soft physics from LHC 
will be available (after year 2005), these might be helpful in reducing some of the uncer- 
tainties in the interaction features and in the shower calculations, as far as the energy 
region of the knee, and above, is concerned. 
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